graphene_matrix_multiply (&s, &u, m);
}
-/* Make a 4x4 matrix that maps
+/* Compute a 4x4 matrix m that maps
* p1 -> q1
* p2 -> q2
* p3 -> q3
* p4 -> q4
+ *
+ * This is not in general possible, because projective
+ * transforms preserve coplanarity. But in the cases we
+ * care about here, both sets of points are always coplanar.
*/
void
perspective_3d (graphene_point3d_t *p1,
#define MAX_ITERATION_COUNT 30
+/* Perform Householder reduction to bidiagonal form
+ *
+ * Input: Matrix A of size nrows x ncols
+ *
+ * Output: Matrices and vectors such that
+ * A = U*Bidiag(diagonal, superdiagonal)*Vt
+ *
+ * All matrices are allocated by the caller
+ *
+ * Sizes:
+ * A, U: nrows x ncols
+ * diagonal, superdiagonal: ncols
+ * V: ncols x ncols
+ */
static void
householder_reduction (double *A,
int nrows,
}
}
+/* Perform Givens reduction
+ *
+ * Input: Matrices such that
+ * A = U*Bidiag(diagonal,superdiagonal)*Vt
+ *
+ * Output: The same, with superdiagonal = 0
+ *
+ * All matrices are allocated by the caller
+ *
+ * Sizes:
+ * U: nrows x ncols
+ * diagonal, superdiagonal: ncols
+ * V: ncols x ncols
+ */
static int
givens_reduction (int nrows,
int ncols,
return 0;
}
+/* Given a singular value decomposition
+ * of an nrows x ncols matrix A = U*Diag(S)*Vt,
+ * sort the values of S by decreasing value,
+ * permuting V to match.
+ */
static void
sort_singular_values (int nrows,
int ncols,
}
}
+/* Compute a singular value decomposition of A,
+ * A = U*Diag(S)*Vt
+ *
+ * All matrices are allocated by the caller
+ *
+ * Sizes:
+ * A, U: nrows x ncols
+ * S: ncols
+ * V: ncols x ncols
+ */
int
singular_value_decomposition (double *A,
int nrows,
return 0;
}
+/*
+ * Given a singular value decomposition of A = U*Diag(S)*Vt,
+ * compute the best approximation x to A*x = B.
+ *
+ * All matrices are allocated by the caller
+ *
+ * Sizes:
+ * U: nrows x ncols
+ * S: ncols
+ * V: ncols x ncols
+ * B, x: ncols
+ */
void
singular_value_decomposition_solve (double *U,
double *S,